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Adviser i Professor  Harry  Bauch 

The  effect  of  small  disturbances  on  laminar  plane 
Couette  flow  Is  studied.  If  the  disturbances  are  Infini- 
tesimal, the  effects  are  described  by  the  Orr-Sommerfeld 
equation.  Chebyshev  polynomials  are  used  to  reduce  the 
problem  to  linear  algebra.  The  resulting  eigenvalue 
problem  is  solved  using  the  generalized  Rayleigh  quotient, 
which  Involves  substantially  less  computer  time  than 
previous  methods.  Accurate  eigenvalues  are  then  computed 
for  higher  values  of  the  parameters  than  has  been  done 
previously.  These  values  further  confirm  the  belief  that 
Couette  flow  Is  stable  under  Infinitesimal  disturbances. 

Using  the  linear  results  as  a starting  point,  the 
effect  of  finite  disturbances  is  studied.  A system  of 
equations  Is  derived  from  the  Navler-Stokes  equations, 
taking  Into  acoount  non-linear  terms.  In  this  case,  the 
flow  becomes  turbulent  for  certain  values  of  the  parameters. 
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1.  INTRODUCTION 


TSiis  paper  is  concerned  wi  til  some  results  about 
the  stability  of  plane  Couette  flow  under  small  distur- 
bances* In  non-dimensional  form,  the  problem  consists 
of  two  dimensional  fluid  flow  between  two  Infinite  planes, 
one  at  y * 1 moving  at  a speed  of  1 to  the  right,  and  one 
at  y ■ -1  moving  at  a speed  of  1 to  the  left*  If  u(x,y,t) 
is  the  horizontal  component  of  the  velocity  of  the  flow, 
and  v(x,y,t)  is  the  vertical  component,  then  u(x,y,t)  ® y 
and  v(x,y,t)  = 0 is  a solution  to  the  problem*  This  is 
the  basic  laminar  flow,  However,  in  practice  laminar 
flow  is  hard  to  maintain,  and  most  flows  tend  to  become 


turbulent  if  disturbed*  Knowledge  of  turbulent  flow  is 
In  a very  incomplete  stage,  with  only  some  general  ideas 
of  the  processes  involved*  Here  we  will  be  concerned  with 
the  somewhat  simpler  problem  of  transition;  with  the 
change  from  laminar  to  turbulent  flow.  This  is  the  question 
of  stability;  that  is,  if  the  flow  is  disturbed,  will  it 
return  to  laminar  flow  (stability)  or  change  to  turbulent 
flow.  Stability  of  course  depends  on  the  amplitude  of 
the  disturbance* 


The  first  major  contribution  to  the  study  of  hydro 
dynamical  stability  can  be  found  in  the  theoretical 
papers  of  Helmholtz,  around  1868.  Rayleigh  and  Reynolds 


2* 

made  further  contributions  near  the  end  of  the  nineteenth 
oentury.  In  particular,  dimensional  analysis  led  Reynolds 
to  the  idea  of  what  is  now  called  the  Reynolds  number, 
a pure  number  inversely  proportional  to  the  viscosity. 
Experimental  evidence  led  Reynolds  to  the  conclusion  that 
stability  breaks  down  when  this  number  exceeds  a critical 
value. 

Explaining  this  behavior  proved  difficult.  The 
key  equation  was  arrived  at  Independently  by  Orr  in  1907 
and  Sommerfeld  in  1908,  now  called  the  Orr-Sommerfeld 
equation.  While  only  valid  for  infinitesimal  perturbations, 
it  gave  a starting  point  for  explaining  transitional 
behavior.  This  equation  was  applied  to  various  types  of 
fluid  flows,  including  plane  Couette  flow,  one  of  the 
simpllest  cases. 

The  Orr-Sommerfled  equation  depends  on  two  parameters* 
the  Reynolds  number  R and  the  wavenumber  o(  of  the  pertur- 
bation. The  first  approach  to  the  equation  was  in  terms 
of  asymtotic  analysis;  in  particular,  picking  o(  fixed 
and  considering  the  limit  as  R approached  infinity. 

Hopf  made  the  first  attempt  in  191^.  which  was  refined 
by  Southwell  and  Chitty  in  1933.  For  Couette  flow,  the 
first  Important  work  was  done  by  Wasow  in  1953  and  Grohne 
in  195^»  who  developed  similar  asymtotic  theories  to  show 
that  the  flow  is  stable  for  fixed  o(  when  R is  sufficiently 
large.  This  has  reoently  been  refined  by  Davey  in  1972, 


! 
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who  also  showed  that  the  flow  Is  stable  for  fixed  H when 
oi  is  sufficiently  large. 

Numerical  solution  of  the  Orr-Sommerf eld  equation 
has  been  difficult.  The  first  solutions  were  obtained  by 
Tollmlen  lnr  1929.  who  also  obtained  a critical  Reynolds 
numbers  Lin  in  19^5  Improved  the  mathematical  procedures 
and  laid  the  foundations  for  a general  expansion  of  the 
stability  analysis.  However,  fitting  the  general  solutions 
to  particular  boundary  solutions  proved  to  be  too  mathe- 
matically complex. 

Numerical  Investigation  of  Couette  flow  had  to  wait 
for  the  development  of  computers.  Grohne  In  195^  developed 
some  results  for  the  speolal  case  oC  ■ 1.  The  first 
extensive  numerical  Investigation  was  carried  out  by 
Gallagher  and  Mercer  in  1962  for  values  of  B ^ 1000, 
which  confirmed  the  hypothesis  that  Couette  flow  Is  stable 
to  Infinitesimal  perturbations.  However,  aso<  S gets 
larger,  the  problem  becomes  more  and  more  singular,  making 
the  numerical  solution  more  difficult. 

Deardorff  In  1963  extended  slightly  the  parameter 
range  covered  by  Gallagher  and  Mercer.  More  recently, 

Davey  has  computed  solutions  for  o(R  up  to  100,000. 

Using  these  results,  Davey  has  detected  a pattern  in  the 
solutions  for  the  range  that  Is  not  covered  by  the  asym- 
totlc  analysis. 

In  this  paper  I will  develop  a more  efficient  method 
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of  computing  the  solutions  of  the  Orr-Sommerfeld  equation, 
giving  aocurate  solutions  relatively  easily  for  large 
values  of  e^R*  The  results  verify  Davey's  analysis* 

The  combination  of  numerical  and  asymtotlc  analysis 
leads  to  the  conclusion  that  Couette  flow  Is  stable  to 
infinitesimal  disturbances  for  all  values  of  the  parameters* 
However,  this  is  not  borne  out  in  practice*  In  a series 
or  experiments  Reichardt  In  1956  was  able  to  maintain 
laminar  flow  only  for  R up  to  about  750*  The  assumption 
le  that  non-linear  effects  cause  the  transition  to  turbulent 
flow* 

Hon-llnear  analysis  Is  relatively  recent*  based 
mainly  on  the  work  of  Keksyn,  Stuart,  and  Watson,  developed 
around  I960*  But  so  far  only  fragmentary  results  have  been 
developed,  and  there  are  essentially  no  such  results  for 
plane  Couette  flow.  In  the  second  part  of  this  paper, 
some  results  of  a non-llnear  analysis  are  discussed*  In 
particular.  If  second  order  terms  are  considered,  a tran- 
sition from  laminar  to  turbulent  flow  occurs.  As  in 
experiments,  the  size  of  the  disturbance  needed  to  produce 
turbulence  depends  upon  the  Reynolds  number. 


2.  DERIVATION  OF  THE  ORR-SOKKERFELD  EQUATION 


> 

1 


i 


The  general  equations  of  a two-dimensional  incom- 
pressible flow  are  the  Navier-Stokes  equations 

(2.1)  ut  + utij.  + vtiy  + px  -(u^.  + ^yy)/R  * 0 

(2.2)  Vt  + UV^  + ^y  ♦ Py  “^XX  + Vyy)/&  « 0 

(2*3)  + vy  = 0 

where  u(x,y,t)  is  the  fluid  velocity  parallel  to  the 
x-axls,  v(x,y,t)  is  tho  fluid  velocity  parallel  to  the 
y-axls,  and  p(x,y,t)  Is  the  kinematic  pressure.  A deri- 
vation is  given  in  Fetchov  and  Criminals  (1967).  For 
plane  Couette  flow,  the  boundary  conditions  are 

(2.4)  u(x,+l,t)  = +1 
v(x,+l,t)  = 0 

The  basic  laminar  solution  is 

(2.5)  u(x,y,t ) « y 
v(x.y.t)  = 0 
p(x,y,t)  **  P(x) 

We  now  consider  small  perturbations  of  the  laminar 
flow,  so  that 

(2.6)  u(x,y,t)  = y + u'(x,y,t) 
v(x,y, t ) « v* (x,y, t) 
p(x,y,t ) » P(x)  ♦ p* (x,y ,t) 

where  u*,  v*.  and  p*  are  small. 
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Substituting  into  the  Navler-Stokes  equations,  we  obtain? 

(2.7)  u*^  ♦ yu*x  + u*u*x  + t*  + T*u*y  + Px  + p*x 

-(u'xx  + U^yy)/B  * 0 

(2.8)  + yr*  + u + p* 

ux  x 7 7 

-(r'xx  + r#yy)/B  “ 0 

(2.9)  u»x  ♦ vry  - 0 

We  ©am  eliminate  the  pressure  term  by  differentiating 
(2.7)  with  respect  to  y,  (2.8)  with  respect  to  x,  and 
subtracting  the  equations.  nils  results  lit 

(2.10)  (u*y  - ▼*x)t  + y(u*Xy  - v*^)  + u*yU*x  ♦ u*u*'Xy 

+▼•  +▼'  U*  •fV’U*'  -U*  V*  -U*V* -V»  ▼* 

yyy  yyxx  xxxy 

»(u*  +U*  — V* • 7*  )/R  ■ 0 

xy  v yyy  u xxy  v xxx  v xyy,/n 
In  the  linear  theory,  we  assume  that  the  perturbations 

are  small  enough  that  we  can  ignore  the  second  order 

terms.  So  (2.10)  becomes 

(2.11)  (u»y  - T*x)t  + y(u*y  - ▼,x)x 

—(u*  + U*  — v*  — V*  )/R  b 0 

1 yyy  u xxy  v xxx  xyy" 

Equation  (2.9)  permits  the  introduction  of  a stream 
function  Y* (x,y,t)  such  that  u*  ■ and  * “Yx* 
Equation  (2,9)  is  then  satisfied  automatically  and  (2.11) 
becomes 


(2.i2)  <ryy  + r**>t + r„  * fxx>x  /J7e*n 

~^Tyyyy  + 2 ^xxyy  + 't' xxxx)/B  “ 0 
We  consider  functions  of  the  form 

^(x.y.t)  ■ j^(y)exp(l^(x  - ct)),  that  is,  periodio 

solutions  in  x and  t,  where  o(  is  the  wavenumber  of  the 
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horizontal  motion  and  o is  an  eigenvalue  to  he  determined. 

Since  we  are  assuming  linearity,  a series  of  solutions  of 
this  type  will  not  Interact  with  each  other* 

Equation  (2*12)  becomes 

(2.13)  - 2e\Z  <f"  + l*R(y  - c)(^'»  -s<2^)  -0 

This  is  the  well  known  Orr-Sommerfeld  equation.  The 

boundary  conditions  (2.4)  become 

(2.14)  j£(+l)  - ^*(+1)  = 0 

In  this  formulation,.  (^(x.y.t)  is  a complex  valued’ 
function.  But  since  the  orlnglnal  Navier-Stokes  equations 
contain  no  complex  quantities,  ^(^(x,  y.t)  + y'(x,y,t)) 

Is  also  a solution,  and  is  real  valued. 

The  problem  has  been  reduced  to  a fourth  order 
ordinary  differential  eigenvalue  equation.  The  imaginary 
part  of  c,  where  c = cr  + ic^,.  becomes  the  relevant  result. 

If  Cj^O,  the  disturbance  decays  with  time.  So  thd  flow  is 
stable  with  respect  to  infinitesimal  disturbances.  If 
0,  the  perturbations  grow  indefinitely  with  time. 

Of  course  once  they  become  of  any  significant  size  the 
linear  theory  no  longer  applies. 


METHOD  OP  SOLUTION 


Various  methods  have  been  used  to  solve  the  Orr 


Sommerfeld  equation,  generally  some  kind  of  finite-diffe- 
rence approximation*  However,  more  efficiency  can  be 
obtained  using  Chebyshev  polynomials* 

Other  sets  of  orthogonal  polynomials  have  also  been 
used  In  stability  problems,  usually  chosen  to  satisfy  the 
boundary  conditions  and  have  some  relation  to  the  Orr- 
Sommerfeld  equation*  The  Chandrasekhar-Reid  functions 
are  the  most  common,  though  other  types  have  been  proposed 
The  advantages  of  the  Chebyshev  polynomials  are  discussed 
In  general  by  Pox  and  Parker  (1972)  and  In  particular  by 
Orszag  (1971),  who  used  them  for  the  similar  problem  of 
plane  Poiseuille  flow* 

Briefly,  Chebyshev  polynomials  are  particularly 
useful  for  obtaining  accurate  solutions.  The  Chebyshev 
polynomial  approximations  are  of  infinite  order,  in  the 
sense  that  errors  deorease  more  rapidly  than  any  power- 
of  1/N  as  N approaches  Infinity,  where  N is  the  number  of 
terms  used  In  the  exranslon*  Hence,  while  other  methods 


may  give  better  answers  for  small  N,  accurate  answers  can 
be  obtained  much  more  rarldly  using  Chebyshev  polynomials 
Another  advantage  is  the  efficiency  with  whioh  the  coeffi- 
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dents  may  be  determined  from  the  differential  equation* 
This  makes  the  process  of  changing  to  a system  of  linear 
equations  much  less  complicated* 

So  we  will  represent  (h  (y)  *=  7L  anTn(y)»  where 
Tff(y)  is  the  n*th  degree  Chebyshev  polynomial  of  the  first 
kind,  defined  by  Tn(cos  •&-)  * cos  ird**  Derivatives  and 
integrals  of  0(y)  can  be  represented  relatively  easily 
by  using  the  relation 

0.1)  cnT*n+1(y)/(n+l)  - dn.2T,n-l(y)/(n“1)  = 2Tn(y) 

where 


r ° 

3tt  ’ 0 


If  n < 0 
If  n = 0 
if  n > 0 


if  n 4 0 
if  n > 0 


The  Orr-Sommerf eld  equation  can  now  be  approximated 
by  a series  of  linear  equations  in-  the  an  by  substituting 
In  the  expansion  for  ^(y)  and  equating  the  coefficients 
of  Tn ( y ) * This  procedure,  used  by  Orszag,  has  certain 
computational  disadvantages,  in  that  the  matrices  produced' 
have  coefficients  of  widely  different  orders  of  magnitude* 
This  Is  due  to  the  infinite  sequence  needed  in  working 
with  derivatives* 

V yV 

For  example,  lffl>»(y)  = £ bT  (y)  = > aT»  (y) 
e*>  n n n 

then  cnbR  ■ 2 pap,  where  «£■*•  Indicates  that  the 

summation  is  In  Jumps  of  two*  However,  the  indefinite 

integral  has  a much  simpler  form.  That  is,  if 


aI 

jPij)  - 1 % /Tn(y)  - £ alT(y) 

/ /t,zO  n'  n W-*  n n , 

represents  the  indefinite  integral  of  (f)  (y) . we  can  find 

V >i  4m  r C/nTti-tlty)  Jm-l  1h>~)  ^ 7 

£an/yy>-4  -rl  ~*+i  " 

~ n / n /rt-d  i 

^/e-/  C/n-j  ~ d/ti-l  fy^' ) 


where  a. 


>«7  s.  0 

0,  so  that 
4/h-l  Cai-I  + l 

2 /VI 


(a^  is  an  arbitrary  constant  of  integration) 

This  is  quicker  to  compute,  and  also  leads  to  a matrix 
whose  coefficients  are  much  closer  in  size* 

Writing  (2.13)  as 

ft  ••••  - 2 4 /- 

i^Rfy^**  -c(2y^-  c^*'  +e(2c^)  - 0 
and  integrating  four  times,  we  obtain 
(3.3)  ^ +<*  *///!& 

-le(lL{{(y0-+zffl0  -o(ZJf/Jyff-  zfjfi  +*Zcfffffi  ) - 0 

Letting //^(y)  =^a2nTn(y) 

///<^<y)  *!Ca3nTn(y) 
and  Jffffl  (y)  «:|>Vn(y) 


we  can  obtain,  as  above 


t,  ^ C/n-i+'A  Gtr\-x+ 


j.  /VI 


.A+A+J+l  Q„-b  +A+A+4 

A*  

^ ao  fan  -z  ) (")-  H~  i~A  +4  ^ 


t+M+^CI^ -i+i+H+kn  Cyi-y-Ui+Al-U^ 

) fa  - + +** H)  (* -61-*” 


and  use  the  relation 


we  can  obtain 


M ~ l 'hti*  •j-A.'+J’  C- Aft- b~{- Q.'f-fl 'f-j* 


^2j;*£*£4£  Qtf-Krt™  +(+fafyC/fi  -tp+im  i fojS.-j-  jUi 

~,  *»*/  fi  &1  ti*1  J^/yi (rY)-2-f'/m')(rr)*  ^ -f-orj  +f)(/v\-l>i-rM+Ji-jL  ) 


(Exception*  If  n * 4,  3=3.  and  m=l=k«h=l, 
the  above  coefficient  is  doubled) 

Substituting  all  of  these  expansions  into  (3*3) » 
equating  coefficients  of  TQ(y),  and  multiplying  by  i, 
we  obtain  the  system  of  eouatlons 
(3.12)  ian  - 21<*2a2n+  io(4ajj  + o(  Rb2n  - 2o<  Ba^ 

- o^Rb1^  - co<Ra21+  o«<2Ra^1*  0 

n'  - 4,5, ... 

Because  of  the  integration,  this  equation  for  the 
cases  n ■ 0, 1,2,3  contains  arbitrary  constants.  Hence 
we  replace  these  equations  by  equations  representing  the 
boundary  conditions.  These  conditions  can' be  simplified 
using  the  relations 

Tn(+1)  = (+l)n 
and1  T»n(±l)  » n2(+l)n_1 

. A' 

thus  0(1)  ■ £ a « 0 

/*>=<?  n 

^(-1)  = £ (-l)na  - 0 
/n-o  n 

0*  (1)  - 2.  n2a_  » 0 

^ m=0  n 

^ «(-l)  - 1l  ^(.ijn-ia  . o 

By  adding  and  subtracting  equations,  we  obtain  the  slightly 
simpler  form 


(3-13)  £ an  = 0 

rr\~0 

#*a  - 0 

h n 

2 

X n a»  ■ 0 

4*-  n2an  . 0 
/*>-/ 

Using  these  four  equations  plus  equation  (3*12) 
for  n«  4,5,...N,  we  obtain  a system  of  N+l  equations 
In  N+l  unknowns.  In  matrix  form  this  can  be  written  as 
XA  - oYA  = 0,  where  X Is  a complex  valued  matrix,  Y Is 
a real  valued  matrix,  and  A » (a^a^  ...aN)  Is  the  vector 
to  be  determined.  So  the  oringlnal  differential  equation 
Is  reduced  to  the  above  algebraic  eigenvalue  problem. 


4.  SOLVING  THE  MATRIX  EQUATION*  THE  LR  ALGORITHM 


There  are  a number  of  standard  methods  for  solving 
an  algebraic  eigenvalue  problem.  Wilkinson  (1965)  contains 
a fairly  complete  description.  However,  the  present 
problem  presents  several  difficulties. 


First,  the  usual  form  for  an  algebraic  eigenvalue 
problem  Is  XA  - cA  =»  0*  An  equation  of  the  form  XA  - cXA  * 0 

Is  normally  rewritten  as  X^XA  - cA  ■ 0.  Unfortunately, 

In  our  problem  X Is  not  Invertible.  The  first  four  rows 
of  X consist  solely  of  ?eros,  since  the  eigenvalue  o does 
not  enter  Into  the  boundary  conditions. 

A method  has  been  devised  by  Gary  and  Helgason  (1970) 
to  get  around  this  difficulty.  A sequence  of  elementary 
column  transformations  is  used  to  put  the  first  four  rows 
of  X into  lower  triangular  form.  This  Is  accomplished 
by  adding  a multiple  of  the  first  column  to  the  succeeding 
columns  so  that  xl9,  x--,,...  all  become  zero.  We  continue 


the  process  for  the  second,  third,  and  fourth  columns, 
permuting  the  columns  If  necessary  to  put  the  largest 
element  in  the  appropriate  column.  Since  the  elements  in 
these  rows  are  known,  this  is  easy  to  program  efficiently. 
The  same  column  transformations  are  performed  on  X. 


The  first  four  rows  of  X remain  zero,  so  If  we  drop  the 
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first  four  rows  and  columns  of  both  X and  Y,  the  new 
system  has  the  same  eigenvalues  as  the  orlnginal  system* 

For  convenience,  this  will  still  be  denoted  XA  - cYA  » 0, 
but  Y is  now  generally  non-singular. 

To  change  this  into  the  form  Y“^XA  - cA  * 0,  programs 
from  Forsythe  and  Moler  (1967)  were  adapted.  The  matrix 
Y is  first  scaled,  so  the  norms  of  the  rows  are  approximately 
the  same  size.  Then,  using  partial  pivoting  to  reduce 
round  off  error,  Y is  decomposed  into  the  form  Y = LTJ, 
where  L is  a lower  triangular  matrix  and  U is  an  upper 
triangular  matrix.  Instead  of  actually  interchanging  the 
rows,  a vector  is  introduced  to  keep  track  of  which  elements 
have  been  used  as  pivots,  saving  computer  time.  Then  for 
each  column  X^  of  X,  the  equation  YZ^  = X^  can  be  easily 
solved,  where  Z ^ is  the  i*th  column  of  Y“^X.  To  insure 
accuracy,  iterative  improvement  is  used.  That  is,  we 
form  a residual  R^  = Xj_  - YZ^  and  solve  the  system 
YZ*i  = R^.  Then  Z + Z*  is  a more  accurate  solution. 

Since  most  of  the  arithmetic  in  this  procedure  Is  in  the 
orlnginal  decomposition,  this  is  inexpensive.  Again,  for 
convenience,  this  new  system  can  be  denoted  XA  - cA  = 0. 

A second  difficulty  is  that  X is  a general  complex 
matrix.  The  most  common  eigenvalue  problem  is  for  a real 
symmetric  matrix,  and  many  of  the  standard  procedures 
do  not  aj?ply  to  complex  matrices.  Wlldinson  discusses 
only  one  method  for  a general  complex  matrix,  the  LR 
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algorithm,  which  Is  relatively  complicated.  To  Implement 
this,  programs  from  Wilkinson  and  Reinsch  (1971)  were 
translated  from  Algol  and  adapted  to  this  particular 
problem. 

To  reduce  round  off  errors,  the  matrix  Is  first 
balanced.  Since  errors  In  eigenvalue  procedures  depend 
on  the  size  of  the  norm  of  the  matrix,  diagonal  similarity 
transformations  are  used  to  reduce  the  norm  of  X.  Air 
Iterative  process  replaces  X by  D“1XD,  until  the  norm  of 
the  i*th  row  and  column  of  X are  of  the  same  order  of 
magnitude. 

Next  the  matrix  is  tranformed  to  upper  Hessenberg 
form,  that  is,  x^j  ■ 0 if  17  J + 1.  While  no  finite 
sequence  of  similarity  transformations  can  reduce  a matrix 
to  upper  triangular  form,  elementary  matrices  can  introduce 
zeros  In  all  places  below  the  subdiagonal  of  the  main 
diagonal. 

Finally,  the  LR  algorithm  is  applied  to  the  upper 
Hessenberg  matrix.  Again,  a sequence  of  similarity  trans- 
formations Is  used  to  reduce  . the  elements  under  the  main 
diagonal.  Since  only  the  elements  In  the  subdiagonal 
are  non-zero,  the  procedure  converges  much  more  rapidly 
than  for  a general  matrix.  As  soon  as  the  term  x2l  becomes 
negligible,  x^  is  an  eigenvalue  for  the  matrix.  The  first 
row  and  column  can  then  be  deleted  and  the  irocedure 
repeated. 
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There  are  some  drawbacks  to  the  above  procedure. 

The  process  Is  long,  being  relatively  expensive  In  computer 
time.  Because  of  the  efficiency  of  the  approximation  by 
Chebyshev  polynomials.  It  does  compare  favorable  with 
other  methods,  but  It  does  compute  a lot  of  Information 
of  no  immediate  value.  In  particular.  In  stability  theory 
we  are  concerned  with  the  eigenvalue  having  the  largest 
imaginary  part,  since  this  will  represent  the  solution 
closest  to  being  unstable.  But  the  LR  algorithm  computes 
all  the  eigenvalues.  Since  in  some  cases  we  will  go  up 
to  IT  = 80,  this  would  be  ?6  eigenvalues  computed  that  we 
have  no  Interest  in. 

So  in  the  next  chapter  we  devise  an  alternate  method 
of  solving  the  eigenvalue  problem,  one  apparently  not 
previously  used  in  stability  problems. 


5.  SOLVING  THE  MATRIX  EQUATION  s THE  RAYLEIGH  QUOTIENT 

The  second  method  of  solution  Is  based  on  a series 
of  papers  by  Ostrowskl  (1958-59)*  These  papers  considered’ 
some  properties  of  the  Rayleigh  quotient* 

As  orlnglnally  defined,  we  assume  that  X Is  a real 

T 

symmetric  matrix,  U Is  a row  vector,  and  U Is  Its  transpose* 
Then  R(U)  = (UXUT)/(UUT)  Is  called  the  Rayleigh  quotient 
corresponding  to  U* 

The  eigenvalues  of  X can  be  characterized  by  the 
extremal  properties  of  this  quotient.  If  c^  Is  the 
largest  eigenvalue  of  X,  then  c^  = max  R(U),  where  this 
maximum  Is  taken  over  all  possible  vectors  U*  If  Is 
the  eigenvector  corresponding  to  c, , then  the  next  largest 
eigenvalue  * max  R(U),  taken  over  all  vectors  U that 
are  orthogonal  to  U^.  The  other  eigenvalues  are  similarly 
characterized* 

In  practice,  what  Is  more  Important  than  the  extremal 
properties  Is  the  fact  that  the  quotient  is  stationary. 

That  is,  Oj,  ■ R(U^)  differs  from  R(U^  + £ W)  only  by 
terms  of  the  second  order  In  £ • So  if  we  have  a vector 
that  is  close  to  the  correct  eigenvector,  the  Rayleigh 
quotient  will  be  quite  close  to  the  corresponding 
eigenvalue. 


18. 


19. 


This  makes  possible  an  Iterative  procedure  for 
determining  any  specific  eigenvalue  c.  We  start  with 
an  approximation  Cq  that  Is  reasonably  close  to  c.  Then 
we  choose  an  arbitrary  vector  Uq.  The  only  requirement 


is  that  Up  cannot  be  orthogonal  to  the  eigenvector  U. 


Then  the  equations 


(5.1)  (X  - ckI)  UTk+1  . UTk 

(5.2)  ck+1  - HktAl’/'Ww' 


are  used  to  get  better  approximations  for  c and  U. 


Ostrowskl  shows  that  the  convergence  Is  cubic,  that  Is, 

(°k+l  “ c)/(ck  “ cP  approaches  a constant. 

If  X Is  not  symmetric,  we  Instead  consider  the 
generalized  Rayleigh  quotient  R(U,V)  = (VXU^)/(VU^) • 

In  this  case,  if  c is  an  eigenvalue,  U Is  the  corresponding 
right  hand  eigenvector  (XU  = cU  ),  and  V Is  the  corre- 
sponding left  hand  eigenvector  (VX  =»  cV),then  R(U,V)  = c. 

The  extremal  properties  of  the  Rayleigh  quotient  no  longer 
hold.  However,  Ostrowskl  shows  that  the  quotient  Is  still 
stationary.  If  the  matrix  X has  only  linear  elementary 
divisors.  This  will  be  true  If  the  eigenvalues  are  distinct, 
which  Is  the  case  in  fluid  mechanics.  Moreover,  the 
convergence  is  still  cubic. 


The  Iterative  process  Is  now  defined  by  the  equations 


20* 


Wx  * V>  - vk 


(5.*)  Ojt+1  - IvAl'/tVl'1  k+1* 


o 


Our  problem  Is  of  the  more  general  form  XU  - cYU  ■ 0* 
The  method  of  Gary  and  Helgason  could  be  used  to  reduce 
this  to  the  standard  form  Y”YXU  - cU  » 0.  However,  the 
Iteration  can  still  be  defined  if  the  equation  is  left  in 
its  orlnginal  form*  The  equations  become 


(5.5)  (x  - ckY)uTk+1  » uTk 


WX  ' \Y>  = Vk 


(5.6)  ck+1  - (Vk+1XUTk+1)/(Vk+1YUTk+1) 


To  actually  carry  out  the  procedure,  we  need  a 
starting  value  cQ.  This  can  be  obtained  by  using  the  LB 
algorithm  for  a smaller  value  of  N,  by  using  an  eigenvalue 
computed  for  nearby  values  of  A and  R,  or  from  the 
asymtotlc  analysis.  In  practice,  the  asymtotic  results 
are  close  enough  for  the  procedure  to  converge* 

The  matrix  X - cQY  is  then  decomposed  Into  the 
product  of  a lower  triangular  matrix  times  an  upper  tri- 
angular matrix.  The  same  decomposition  can  be  used  in 
obtaining  both  U and  V,  since  the  equation  (5*5)  becomes 


(5.7) 


T 

(Lower) ( Upper )U  k+1  « 


k+1 


3 


(Upper  )T  ( Lower ) TVT. 
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where  (Upper)T  is  lower  triangular  and  (Lower)T  Is  upper 
triangular.  Since  frost  of  the  computer  time  In  solving 
a system  of  equations  is  required  for  the  decomposition, 
this  procedure  results  In  a noticeable  time  savings. 

To  help  Insure  that  UQ  and  V0  are  of  the  proper  form, 
they  are  defined  implicitly  by  making  (Lower )“^UT0  and 
((Upper)T)“1VT0  both  equal  to  the  vector  whose  components 
are  all  ones.  The  vectors  then  depend  on  ths  matrix,  and 
the  values  computed  for  and  are  unlikely  to  be 
orthogonal  to  U and  V. 

For  the  first  several  iterations,  only  the  equations 
(5»5)  &re  used.  This  allows  the  eigenvectors  to  converge 
relatively  closely  to  the  proper  values.  Again,  this 
involves  only  the  oringlnal  lower  upper  decomposition, 
so  the  procedure  is  efficient.  Then  the  full  iteration 
can  be  used  until  the  eigenvalue  converges.  Generally, 
the  eigenvalue  will  be  accurate  to  four  decimal  places 
after  only  three  or  four  additional  Iterations. 

The  two  methods  were  compared,  and  they  do  lead  to 
the  same  numerical  answers.  However,  the  generalized 
Rayleigh  quotient  is  much  more  efficient,  particularly  for 
large  N.  For  example,  if  N * 70,  the  generalized  Rayleigh 
quotient  was  more  than  five  times  as  fast  as  the  LR 
algorithm. 


6.  NUMERICAL  RESULTS 


All  existing  evidence  tends  to  show  that  plane 


Couette  flow  is  stable  under  infinitesimal  disturbances 


However,  because  the  eigenvalue  c is  a function  of  two 
parameters,  the  wave  number  o(  and  the  Reynolds  number  R, 
it  is  not  easy  to  systematically  cover  the  entire  •(  ,R 
plane.  Moreover,  as  the  product  of  and  R becomes 
larger,  the  Orr-Sommerf eld  equation:  becomes  more  and  more 
singular,  and  hence  harder  to  solve  numerically.  And 
any  instability  would  be  expected  to  occur  for  large  values 


The  most  recent  results  on  plane  Couette  flow  are 
due  to  Davey  (1973).  who  used  a combination  of  asymptotio 
and  numerical  methods.  Based  on  the  analysis  of  Wasow 
and  Grohne,  Davey  studied  asymptotically  the  most  familiar 
oase,  where  >TfT  is  much  larger  than  , and  is  much 

larger  than  1.  Using  Airy  functions,  he  arrives  at  air 
estimate  for  the  eigenvalue  c * cp  + icA  having  the  largest 
Imaginary  part,  given  by 

//» 

(6.1)  or  » 1 - 4.1287/(<*R) 

o>  « -<*/R  - 1.0625/(«<R)//7 

-Yt 

He  estimates  the  error  as  being  of  order  R • 

This  indicates  that  for  a given  fixed  , as  R 


approaches  Infinity,  the  real  part  of  the  eigenvalue 
approaches  one  and  the  Imaginary  part  approaches  zero* 
However,  the  Imaginary  part  Is  always  negative.  Indicating 
that  the  flow  Is  stable* 

Table  1 shows  the  relative  error  In  the  approximation 
for  = 1 and  different  values  of  R.  As  R gets  larger, 
the  asymptotic  results  get  closer  and  closer  to  the  actual 
values  * 

TABLE  1 

COMPARISON  OF  ASYMPTOTIC  AND  COMPUTED  VALUES  FOR  «<  = 1 


R 

r— 1 — 

Asymptotic  Results 

Computed  Results 

■ 

Error 

500 

.4800  - .13581 

.5092  - .15451 

12.24* 

1000 

.5871  - .10731 

.6053  - .11921 

3.52* 

5000 

.7586  - .06231 

•7646  - .06661 

.96* 

10000 

.8084  - .04921 

.8122  - .05211 

.59* 

50000 

.8879  - .02881 

.8892  - .02981 

00 

H 

• 

100000 

.9110  - .02291 

.9119  - .02361 

.12* 

200000 

.9294  - .01021 

.9299  - .01861 

.07* 

For  the  larger  values  of  R,  the  Orr-Sommerf eld 
equation  Is  more  difficult  to  solve.  But  using  Chebyshev 
polynomials,  this  simply  means  that  N,  the  degree  of  the 
largest  Chebyshev  polynomial  used,  must  be  Increased. 
Thus,  If  R a 500,  N = 40  gives  sufficient  accuracy,  while 
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if  R » 200000,  we  must  choose  N - 80.  If  larger  values  of 
H or  more  than  four  place  accuracy  is  required,  the  value 
of  N then  must  be  chosen  larger.  A larger  value  of  N 
does  require  more  computer  time,  but  the  process  is  still 
quite  efficient. 

The  asymptotic  results  can  be  Improved  by  noting 
that  the  answers  given  are  always  numerically  smaller 
than  the  actual  answers.  Using  Davey's  error  estimate, 
we  oan  create  a new  approximation 

(6.2)  o*l-  4.l287/(«<R)/^+  R~ 

~X/J 

o , - -°</R  - 1.0625/(«<R)  - R 

• IZ 

As  R approaches  infinity,  R ^approaches  zero,  so 
this  doesn't  effect  the  asymptotic  theory.  However,  it 
does  give  substantially  better  numerical  results  for 
smaller  values  of  R,  as  can  be  seen  in  Table  2,  where  the 
relative  error  is  recomputed  using  this  new  approximation. 
Using  these  modified  results,  the  generalized  Rayleigh 
quotient  normally  converges  rapidly  to  an  answer. 

Davey  also  studied  the  less  familiar  case  when 
R is  large  and  o(  is  much  larger  then  /~lT.  Again 
using  Airy  functions,  he  arrives  at  an  estimate 

(6.3)  cr  - 1 - 2.0249/H  *)//j 

//* 

c.  ■ — °VR  - l.l69l/(°<  R) 

~r/> 

where  the  error  is  of  order  «*(  • This  indicates  that  for 

a given  fixed  R,  as  c<,  approaches  infinity,  the  flow 


remains  stable. 


25. 


' 


TABLE  2 

COMPARISON  OP  THE  MODIFIED  ASYMPTOTIC  AND  COMPUTED  VALUES 

FOR  = 1 


R Asymptotic  Results'!  Computed  Results  Error 


1000 

5000 

10000 

50000 

100000 

200000 


.5092  - .15451  ' 4.80# 

.6053  - .11921  1.36# 
.7646  - .06661  .36# 
.8122  - .05211  .22# 
.8892  - .02981  .08# 
.9119  - .02361  ;o5# 
.9299  - .01861  .0^# 


The  asymptotic  theory  Is  only  valid  when  either 
or  H Is  large.  However,  If  both  parameters  are  small, 
the  Orr-Sommerfeld  equation  is  relatively  easy  to  solve. 
Gallagher  and  Mercer  (1961)  examined  numerically  the 
eigenvalues  for  <*R£-i000.  For  these  small  values,  the 
flow  Is  stable. 

The  other  case  that  Is  not  covered  by  the  asymptotic 
theory  is  when  oC  and  /h  are  of  the  same  order  of  magnitude 
This  means  that  the  terms  of  the  Orr-Sommerfeld  equation 
are  of  the  same  order  of  magnitude,  and  no  asymptotic 
simplification  Is  possible. 
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To  explore  these  values,  Davey  introduces  new 
variables 
( 6A ) X =c^H'Vi 

r - -r^>1 

Then  the  second  part  of  (6.1)  becomes 

(6.5)  Y *=  X + 1.0625X  + 0(X7') 

Where  X is  small  (that  is,  o(  is  much  smaller  than/ft) 
and  R is  large.  The  second  part  of  (6. 3)  becomes 

(6.6)  Y e X 4*  1.1691X  + 0(X  *) 

where  X is  large  (that  is,  is  muoh  larger  than  /r) 
and  R is  large.  This  suggests  that  for  large  values  of 
R,  when  is  not  too  small, 

(6.7)  Y = X + HU) 

where  P(X)  is  a slowly  varying  function  of  X whose  value 
remains  close  to  I.O625  and  1.1691*  If  this  is  true, 
then  for  large  enough  R,  Y is  Just  a function  of  X.  If 
Y is  always  positive,  will  be  negative,  and  the  flow 
will  be  stable  to  Infinitesimal  disturbances  for  all 
values  of  and  R. 

Davey  claims  that  this  is  the  case  for  R_>200. 

While  not  giving  the  exact  numerical  data,  he  states  he 
has  computed  the  eigenvalues  for  values  of  «<  R up  to 
100000,  and  the  values  for  Y depended  solely  on  X.  The 
resulting  curve  is  given  in  Figure  1,  which  also  shows 
the  curves  predicted  by  the  asymptotic  equations  (6.5) 
and  (6.7)  for  comparison. 


Pig.  1.— Graph  of  Y versus  X 


^ yzx+i.uHix"'* 


/;  X+I.OHS  X 


To  check  this  result,  eigenvalues  were  computed  for 
values  of  oC.  R up  to  250000.  The  results  for  R = 200, 

R = 1000,  R = 3000,  and  R = 5000  are  given  In  the  following 
tables.  The  results  show  a high  degree  of  consistency. 

For  the  larger  values  of  °<  R the  results  do  not  match 
exactly.  However,  In  this  case  a small  error  In  cA  can 
lead  to  a fairly  large  error  In  Y.  And  the  general 
principle  that  Y depends  on  X seems  to  be  well  established. 
Any  variance  In  the  values  of  Y appears  to  be  quite  small 
compared  to  the  size  of  Y. 
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TABLE  5 


VALUES  OF  X,  C.  AND  Y FOR  R * 3000  AND  N * 70 


.9296  - .03331 


28.28  .9386  - .03191 
35.36  .9448  - .03101 
42.43  .9491  - .03061 


2.529 

2.356 


2.167 


.9527  - .03101 


2.189 





To  check  this,  some  additional  values  were  run  for 


4 , then 


For  this  case,  if  R 


Y « 2.25^,  which  exactly  matches  the  value  for  R = 200  and 


7,  then  Y = 2.l8l,  which  is 


much  closer  to  the  proper  value.  And  if  R = 5000  and 
X = .4,  then  Y = 2.258,  which  is  slightly  closer  to  the 


expected  value.  So  if  N is  chosen  large  enough,  the 
values  of  Y do  appear  to  agree. 


For  a given  fixed  value  of  R,  the  minimum  value 


of  *c.  occurs  between  X 


c.  will  approach  zero,  but  it  will  always  be  negative 


The  belief  that  Couette  flow  is  stable  under  infinitesimal 


disturbances  is  further  confirmed 


; 

. ( 


f 


7.  THE  NON-LINEAR  ANALYSIS 

The  linear  case  predicts  stability  for  all  values  of 
and  B.  However,  as  has  been  mentioned,  Reichardt  (1956) 
was  able  to  maintain  laminar  flow  only  for  values  of  R up 
to  about  750* 

This  could  be  due  to  various  causes.  First,  the 
apparatus  could  only  approximate  an  infinite  channel.  Also, 
the  flow  in  an  experiment  is  three  dimensional,  and  vortices 
perpendicular  to  the  direction  of  the  flow  could  affect  the 
results.  Finally,  the  perturbations  are  of  course  finite, 
and  there  is  no  exact  way  of  measuring  their  amplitude  or 
wave  number.  In  this  section  of  the  paper  I will  try  to 
take  into  account  this  last  effect,  by  assuming  a non- 
linear solution. 

The  main  approach  to  this  type  of  problem  was  deve- 
loped by  J.  T.  Stuart  and  J.  Watson  in  a series  of  papers 
in  the  early  1960*s.  We  will  use  here  In  general  the  notation 
of  Reynolds  and  Potter  (1967),  vrho  modified  somewhat  the 
approach  of  Stuart  and  Watson. 

We  will  again  consider  the  Navler-Stokes  equations 

(7.1)  ut  + uux  + vuy  + px  - (uxx  + Uyy)/R  = 0 

(7.2)  Vt  + UVX  + Wy  + Py  - (V^  + Vyy)/R  = 0 

(7.3)  ux  + vy  = 0 
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with  the  boundary  conditions 

(7.4)  u(x,  +1,  t)  = +1 
v(x,  +1,  t)  * 0 

Then  we  attempt  to  solve  this  in  terms  of  the  basic  linear 
stream  function  '//(x,y,t)  = ^(y)  exp(i<^(x  - ct))  and 
its  harmonics,  under  the  assumption  that  the  linear  solution 
is  a reasonable  approximation  to  the  actual  solution. 

In  the  general  case,  the  amplitude  of  the  distur- 
bance becomes  important,  so  we  will  try  to  isolate  this 
term.  The  linear  solution  can  be  rewritten  as 

(7.5)  ^(x.y.t)  = ft  (y)  exp(iK(x  - crt))  exp^cjt) 

In  this  case,  the  magnitude  of  the  disturbance  is  determined 
by  <f>  (y)  and  exp^c^t),  which  are  Independent.  The  Orr- 
Sommerfeld  equation  determines  (p  ( y)  only  up  to  an  arbitrary* 
multiple,  so  it  can  be  normalized  in  any  convenient  fashion. 
This  determines  the  relative  size  of  the  disturbance  for 
different  values  of  y,  but  the  absolute  size  of  the  distur- 
bance is  immaterial.  The  term  exp(°\c^t)  determines  how 
the  disturbance  behaves  with  respect  to  time,  that  is, 
whether  it  will  die  out  or  increase  to  the  stage  where  the 
linear  theory  no  longer  applies.  As  far  as  we  are  concerned, 
this  is  the  important  part  of  the  solution. 

To  try  to  Isolate  a similar  function  in  the  general  case, 
we  introduce  a change  of  variables 

(7.6)  x + lv  t 
A - A(t ) 

h W rn  \JJ  (p) 
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Here  -0-  represents  the  periodic  part  of  the  solution,  °C 
as  usual  represents  the  wavenumber,  and  U>  represents  the 
frequency  of  the  solution.  In  the  linear  case, 

O’  - x - cK  opt  and  U> = -o(cr.  In  the  non-linear  case, 
the  frequency  UJ  will  depend  on  the  size  of  the  disturbance. 
So  Ut  Is  viewed  as  a function  of  a yet  to  be  defined 
amplitude  function  A,  which  depends  only  on  the  time  t. 

By  analogy,  the  function  A(t)  will  show  how  the  disturbance 
changes  for  any  fixed  y with  respect  to  time.  In  the  linear 
case,  A(t)  = expfc^c^t). 

With  this  change  of  variables,  the  Navier-Stokes 
equations  become,  with  u and  v considered  as  functions  of 

£ + rf*. 

(7.6)  ^ + [ w + 1 ^ 


(7.9) 


+ - a C** 

H-  ^ -I^T  - ^ 


r , i 1-  ^ 

L * J-O 


As  before,  we  eliminate  the  pressure  term  by  differentiating- 
(7,7)  with  respect  to  y and  multiplying  (7,8)  by 
ck  and  differentiating  It  with  respect  to-0-,  and  sub- 
tracting the  equations.  This  results  in 


p i. 

(7*io)  0 MPf  +E  w A ^ 3^  /yP*- 

+ «%rii  -t-ff  ij.  + rjy-  + ‘*'j£dr 

~r[°<  -f-^-j_^dA  J>*r 

^ *y  Fpj  *3*  jAj-fi 

-*["++ ft# 

+«$7*J=  a 

Again  we  Introduce  a stream  function  /'(A.y,^-) 


defined  by 


<7.n)  |p  = ^ ~ r 

Also,  for  convenience  we  Introduce  the  notation 

(7.12)  -f  (A,  P,  ■&■)  - Yjfz  + * * 7+*- 

Substituting  and  multiplying  by  C*s  , equation  (7»10)  becomes 

(7.U)  +t!M  +$£]§£ 


fy- 


M r ilf 


The  continuity  equation  (7.3)  becomes 

(7  *,  All  - - 0 

<7  ’ d+»r  ~ u 

so  It  is  automatically  satisfied. 

The  boundary  conditions  (7**0  become 

(7-15)  §£(a>  ± !>  *'>~0 
|p  <'4,H  «-)  - ± * 


35 


We  can  now  expand  the  stream  function  as  a Fourier 


series 


<7.i6>  H(A,y,+)  - X [ 

^CH)  /A) 

where  i'  yis  the  complex  conjugate  of  7 • So  while 

the  fma.v)  are  complex  valued,  ^ (A,y,  ^")  Is  real  valued* 

To  reduce  somewhat  the  mathematical  complexity  of  the  problem, 

we  will  terminate  the  expansion  at  k » 3*  Since  the  linear 

case  has  been  chosen  as  satisfying  the  Orr-Sommerfeld 

equation,  hopefully  this  will  already  be  reasonably  accurate* 

Note  that  ^(A.y)  can  be  chosen  real  valued,  since  any 

Imaginary  part  will  cancel  out. 

So  we  have  • 

(7.1?) 


-h 


We  introduce  the  notation 


This  allows  us  to  write 


o + , ?ro\ 


(7.19)  = 2^-’*+ 


We  substitute  (7.17)  and  (7*19)  Into  (7*13)  and  equate 

, i ^ 

the  coefficients  of  the  constant,  £-  , A*  , and 

terms.  We  will  Ignore  all  terms  of  order  greater  than  three. 


I 
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The  last  of  these  four  equations  will  not  be  required 
in  our  procedure,  so  we  will  Just  give  the  first  three 


equations  yr/)  _ 

O)  &4  & - o 2 ~h  " b'y 

(7.20)  TA  C /»  . 

-it  = a 

<7-21)  1&  $+*[”  #1*”+ 1 * 7]T 

+ tf£2t0-y*  rfa4^ 

+ i P>  $f- * *»  0'L* **  Is  * 

»*>&$?+["+*#& ]**  *M  s,„ 


The  boundary  conditions  (7.15)  become 
1 \ list*)  - 


(7.23)  9.  fAy^O-  ^ 

yfl){4/  ±/)  r {A,  * 0-f 

The  fact  that  the  equations  are  non-linear  and  coupled 
makes  the  solution  difficult.  However,  assuming  that  the 
amplitude  is  small,  we  can  seek  a solution  as  a power  series 
in  A.  We  require  that  the  solution  for  infinitesimal 
amplitude  reduce  to  the  Orr-Sommerfeld  wave,  and  the  solution 
for  zero  amplitude  to  the  basic  laminar  solution.  This 
suggests  seeking  a solution  in  the  form 


37. 


ru a,?)--  i a"* 


* A> 


(7.24) 


/T) 


A 


We  will  show  that  this  satisfies  the  above  conditions. 
Again,  we  will  truncate  the  expansion  at  n * 3,  and  k 
will  only  have  the  values  0,1,2,  and  3. 

We  also  expand  as  a power  series 


d A.  - & 

(7.25)  M 


to) 


A -h  otl>A''+  aMA'-+- 


For  Infinitesimal  A,  the  amplitude  behaves  exponentially, 
as  in  the  linear  theory. 

Finally  we  represent 

+ ditclA-L  ">+  -t  / ‘*A V — 

* J A ~TT  ~ 'V  1 


(7.26)  yf  Zf 


To  simplify  the  equations  somewhat,  we  Introduce 
the  notation 


(7.27,  z'*' a?- «'£*)  <p ' 


This  allows  us  to  write 

'*>  s t A*2r*;0,> 


(7.28) 


/n-  4 


We  substitute  these  expansions  into  the  equations 
(7.20),  (7.21),  and  (7.22),  and  equate  the  coefficients 
of  the  powers  of  A.  We  will  go  only  up  to  the  A^  terms, 
Equation  (7.20)  becomes  the  four  equations 


■: 


(7.29) 


3 8 

(7.30,  7 4MZ('',)-  0 

+1*  **!”+  z 

+i  0$,t',)?tr,l-*lf/l'nl>*r'll)+'<'  <P"',)D2"',>-ks>  Z10'^-  6 

,7.32)  i a'»2"'W>^+2«'»^l>- 

of(>i»  « 0/,,/  ' 1 

_ u ^ b £ d,l)+vy>'»t>  * 

Equation  (7.21)  becomes  the  three  equations 
(7.33) 

.2,1 Z ° = 0 


(7.34)  2 Q ‘ i ‘ -h  Q 


n;o 


, b -£W  i<r,n  + 2i  />  4>“’'‘)irr^'‘-b2»  b<tte,''}2 

(7.35)  2 Qr0 Sr'1'1*1  + Ci'*1  £ f,rn 

, i +* *■"> 2 ">» + i ^ rr'° 

DfWi*”-  * o <t>  te',,,> 

-ib  <pai*> bi,,i*)+i$0Z)t)i,i'i-li'4a:'>>l>%0''>-  kfi-"1)*  ' _<5 


Equation  (7.22)  becomes  two  equations,  but  only  the  first 


* 

I 
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(7.36  0***2 

+,t><t>"''''2"‘')-A,4>n:'>b  1»  <p'^i>  il * - o 

will  be  needed. 

The  boundary  conditions  (7.23)  become 


(7.37) 


k = 1,2,3 

n *s  k • • «3 


/>  (±l)~  0 >r-  1.2.3 

We  do  not  have  boundary  conditions  ^^"and1 

$ ' • This  depends  on  the  normalization.  What  will  be 
most  convenient  is  to  set  the  bulk  average  velocity 
equal  to  zero,  that  Is,  u(x,y,t)dy  * 0.  This  Is  true 
for  the  laminar  case.  Then 

±l;  %dp--o 

rtA,i,  #■)- 1 rA,-i,-e-)-o 

- ro)fA,  -/)-  0 


This  can  most  easily  be  satisfied  by  choosing 

(7.38)  (p  10 ' ^ (±0=  0 


n = 0,1, 2, 3 


We  will  now  try  to  solve  out  set  of  equations. 
Equation  (7.29)  with  its  boundary  conditions  can  be 


written  as 


40 


(709)  D * f y-)=  o 


(% D-o  b f± /)  - ± "T 


This  has  the  solution 

(7.40)  = TJ  ? ~ If 

Note  that  If  the  amplitude  of  the  disturbance  A ■ 0, 
we  get  V'  (A,y,-£-)  . 2^(y)  « * y2/2  - <*/2  and 
u(x,y,t)  * y and  y(xfy,t)  - 0,  This  Is  the  basic  laminar 
solution,  as  desired*. 

Equation  (7.30)  becomes 

(701)  = 0 

tf,/’;* (to  - b <P n (±<):o 

So  the  solution  Is 

(702)  (y-b-O 

Equation  (7.33)  becomes 

(7.43 )(£>*'- 0 

d>w>&,)  - b<p",» /±,)Z $ 

This  Is  the  Orr-Sommerfeld  equation  with  o - (-b^V  la'*) /<<  , 


that  Is 


(7.44) 


Jr'»  - - * 


C/ 


CO)  = u C v 


If  we  assume  an  Infinitesimal  disturbance,  and 

Ignore  all  but  the  terms  that  are  linear  in  A,  then 
ft) 

dA/dt  = a A a cK.  c^A  or  A = expH  c^t ) , Also* 

= °(  y2/^  - °(  /4  and  exp(o<  c^t ) • b ^ 

will  later  be  shown  to  be  zero,  so  if  we  assume  that  ^ 
is  a constant  (that  Is,  the  frequency  of  the  wave  doesn't 
change),  then  ^ = b^  = -°(cr.  This  gives  us 

t(A,;,+)  '-$?■ - i + 

So  our  formalism  Is  an  extension  of  what  has  been  done 
in  the  linear  case. 


must  be  normalized  in  some  fashion.  We 


will  use  the  standard  normalization 

(7*45)  <p°/n  (O)-} 

This  will  implicitly  determine  the  function  A(t),  which 
is  the  other  part  of  the  solution  relating  to  the  size 
of  the  disturbance. 

Continuing  our  discussion  of  the  separate  equations, 
equation  (701)  becomes 

(7.46)  X »*  + 

i r/iAa-s)  i t";o- ; $ "''U i # 

$ (x  /;--/>  0 ,o;  {z/)-o 

This  is  an  inhomogeneous  differential  equation.  If 
& ' is  known,  this  can  be  solved  for  0^'*} 


Similarly,  Equation  (7*36)  becomes 


(7.47 ){(£+#-*><*(?- 
- „ib<P"'')(bx-*i)0r,'l)-  * 

0'x,*>  f±/)  - 0-6 


This  is  also  an  inhomogeneous  differential  equation. 

The  solution  of  the  last  three  equations  (7*32), 

(7*3*0 » and  (7*35)  requires  some  knowledge  of  the  constants 

a^#  a^5,  b^*,  and  b^.  For  example,  ( 7 • 3^ > becomes 

- . ..j  . „ / . -lr',OXoal‘OrtX^')i0n'O 

(7.48)£/7> -•< ) ' ¥ ~ ^ '‘6  <> 

_ J>  ( a,°0  i 6rm ) (f-**-)4ni0 

(±0-  / i))- 6 

The  difficulty  is  that  there  is  not  enough  information: 
to  specify  the  values  of  a^and  b^.  So  we  must  make 

some  additional  assumption  to  guarantee  a unique  answer. 

(d) 

This  would  not  be  a problem  is  a = °\  c^  « 0.  In 
this  case  the  operator  acting  on  ^ would  be  the  Orr- 
Sommerfeld  operator.  So  the  adjoint  homogeneous  equation 
would  have  a solution,  and  (7.48)  could  not  be  solved 
unless  its  right  side  was  orthogonal  to  this  solution. 

This  consistency  condition  would  determine  a^and  b^. 

This  can  be  done  in  the  related  case  of  plane 


Polseuille  flow,  where  this  method  has  been  previously 
applied.  Here  there  is  a stability  curve  for  the  linear 


I 
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case  in  the°K  -R  plane,  where  = 0.  On  this  curve, 
the  constants  are  determined  by  the  orthogonality  condition* 

Near  the  neutral  curve  the  same  method  is  used  based  on 
continuity.  In  this  case  a^  turns  out  to  be  zero,  and 
dA/dt  = a^  + ....  So  the  sign  of  a!  ^ determines 
whether  or  not  the  flow  is  stable  and  the  disturbance 
dies  out  (a  < 0)  or  unstable  and  the  disturbance 
grows  ( a (i)  > 0 ) * 

In  the  case  of  Couette  flow  the  Imaginary  part  of 
the  eigenvalue  is  never  equal  to  zero.  So  we  will  apply 
an  alternate  method,  also  suggested  by  Reynolds  and  Potter, 

Instead  of  insisting  that  the  linear  flow  be  at  the 
transition  point,  we  will  insist  that  the  non-linear  flow 
be  at  the  transition  point,  that  is,  dA/dt  = 0;  So  A(t) 
will  be  a constant,  and  the  disturbance  will  neither 
grow  nor  decay.  The  solution  will  be  steady  state  with 
frequency  U)  , 

To  solve  this  problem,  we  expand 

(7.49)  LO  - U>  (o)  i-  A M l°  f A ^ ' ' ' 

While  is  a real  number,  we  will  view  it  as  complex 
as  a mathematical  convenience.  The  physical  problem 
will  only  have  solutions  for  certain  values  of  A. 

Our  procedure  is  the  same  up  to  the  equations 
(7.20),  (7.21),  and  (7.22).  We  now  set  dA/dt  equal  to 
zero  and  substitute  in  the  expansion  for  , This  creates 
minor  changes  in  the  set  of  equations  (7.29)  thru  (7.36). 


mmmr*'’ 


pm 


p»'_j 
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Equation  (7.29)  remains  the  same,  so  the  solution 
for  In  (7.40)  also  remains  unchanged. 

The  first  term  of  equation  (7.30)  disappears, 
but  the  solution  for  ^^given  by  (7.42)  Is  still  the 
same. 

Por  equation  (7»33)  the  term  &/o) z (>/ ^ drops  out, 
and  Instead  of  the  form  given  by  (7*43)  we  obtain  the 
similar  form 


(7.50 *>)}  <P",n - O 

r /)  (f>  (!?>($.  ()  - 6 


This  is  still  the  Orr-Sommerf eld  equation,  t t now 

-°(o,  where  c Is  the  eigenvalue  of  the  Orr-Sommerf  eld 
equation. 

o>  — , -//) 


For  equation  (7.31),  the  af  ' and  a‘/;  terms  drop 
out,  and  (7.46)  simplifies  slightly  to  become 

ft”'0 


~h  v 


<p^'x)  C±/)  - b0r^]  C±/)-6 


(*> 


Similarly,  the  a terms  drops  out  of  (7.36),  and 
(7*47)  becomes 

> to)  v , ^ „ s \ ~y  j.  ( 


<7.52 )£(£-+*'')  R (?+  5^- 


I 


The  left  Bide  of  this  equation  Is  the  Orr-Sommerf eld 

operator , but  with  replaced  by  2<^  . 

Equation  (7.3*0  was  the  first  of  the  equations  that 

to) 

presented  difficulties.  But  now  with  the  a term  removed. 


(7.^8)  can  be  rewritten  as 


- i 0tr,n 

(£  l):b<P"'*(  i lY-O 


)?  4 


//;j) 


The  left  side  of  this  equation  is  Just  the  Orr-Sommerfeld 

skD'.l) 

operator  acting  on  Y . We  know  that  the  corresponding 


homogeneous  equation  has  a non-trivlal  solution 


0n'i'\ 


So  the  adjoint  homogeneous  system,  which  has  the  form 

(7.54 4>:  ° 

- bpr±/)~£> 


has  a non-trlvlal  solution.  We  will  call  this  solution  v(y). 

The  Inhomogeneous  system  will  have  a solution  only 
If  the  consistency  condition 


(7.55)  4/ 


t ¥r',)(y)']  Vty)dy  = 6 


is  satisfied.  This  Is  unlikely  to  be  true  unless  U>  ' ■=.  o. 
Then  ^^can  be  any  multiple  of  For  convenience 


m 
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we  will  choose  the  solution 
(7.56)  O 


This  leaves  two  more  out  of  our  orlnglnal  set  of 
eight  equations  to  consider.  Equation  (7*32)  now  reduces  to 


(7.57) 


t'W  (y)-0 


(f>  f0> 3>  (±  /)-  D 


so  the  solution  Is  simply 

(7.58)  ( y)  ~ 0 

The  last  equation  (7*35)  can  be  written 

(7.59)  *4r  ) ( if- ^ $ ' 

(p<'f/v  (±/)z  1>  (pn,3> 

The  left  side  Is  still  the  Orr-Sommerfeld  operator,  so 
we  have  a consistency  condition 

(7.60)  [ {£>*-  4^)  IT/f)  dy*  - 

J ‘ vty )[-iDd> 10 ' 5 4 *)  4>f>/  - ib4>0,,)({ 


(!> 2-  <t*)  #tr'n  + 2 b ’ . 

_ 0e,r}#/~ 


fX  ) 

This  determines  the  value  for  Uf  , which  is  what 
we  are  actually  after.  Additional  values  can  be  gotten 
using  the  same  procedures,  but  will  turn  out  to  be 

[U\ 

zero,  and  to  obtain  would  require  much  more  numerical 


work, 


So  we  have  the  relation 


(7.61)  U>  - lMCDl  ~j-  WC^A 

Since  W and  A are  both  real,  this  means 

(7.62)  U>  : Re ~t  Rea)  (W{  ) A 

0 - 0*>w)  + J ”a3  ^ 

This  determines  the  size  of  the  amplitude  A,  which  is 


given  by 


o<  Cl 


(7.63)  fi=  i z„af(£7») 

Since  c^  is  always  negative,  Imag(tf/2?)  must  be 
negative  for  A to  be  defined.  If  Imag(u/^)  is  positive, 
there  is  no  steady  state  solution.  The  flow  is  either 
stable  with  respect  to  any  second  order  disturbance,  or 
else  any  such  disturbance  leads  to  turbulence. 

If  A is  dofined,  it  gives  the  transition  level 
below  which  the  disturbance  is  damped  out  (infinitesimal 
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case)  and  above  which  it  arrows  (instability). 

Some  idea  of  the  physical  significance  of  A can 
be  seen  by  considering  the  center  line  fluctuations. 
We  know  that 


(7.64)  v(A,y.-*)  - ~ fe-  ~ ~u^  ^ 'tvY  X 


Since  A is  small,  as  a first  approximation  we  can  Ignore  any 
power  of  A larger  than  the  first.  Then  (7.64)  becomes 

-l  A 


-<&• 


(7.65)  v( A.y,-^)  « 


or 


(7.66)  y(A,y,  a-)  . 2 A 

Since  we  choose  a normalisation  j^^(0)  « 1,  we  get 


(7.67)  v(A,a-£0  * I A & 


So  A is.  approximately  half  the  amplitude  of  the  vertical 
center  line  disturbance. 


i 

i 


8,  METHOD  OF  SOLUTION 

As  In  the  linear  case,  we  approximate  the  various 
functions  by  a sunt  of  Chebyshev  polynomials.  The  linear 
eigenvalue  c = and  the  linear  eigenfunction 

are  computed  using  the  generalized  Rayleigh  quotient 
iteration. 

To  find  the  adjoint  eigenfunction  v.  we  use  the 
same  Inverse  iteration,  except  the  eigenvalue  is  now  kept 
fixed.  Equation  (7*54)  Is  Integrated  four  times  to  obtain 

(8.1)  j,\r  ~ ib*.* ff\r  * ffj f 1 r -f  u R //  y*  ^ 

- R Wf  y \r-  Rc  ff\r  fflfir-d 

N 

Letting  v(y)  « a^fy)  and  substituting  into 

(8.1)  we  get  the  system  of  equations 

(8.2)  -b 

- R - * % C Om  ~h  ^ R CCfl,  = 0 

* 4,5. » 

The  only  difference  between  this  system  and  the  system 
for  the  Orr-Sommerf eld  equation  (3.12)  is  that  the  term 
-2p(  Ra  n is  missing.  After  making  this  change,  v(y) 
can  be  computed  using  the  same  procedure  as  for 
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To  solve  the  other  equations,  we  need  subroutines 
for  talcing  Integrals  and  derivatives  and  for  multiplying 
The  formula  for  Integrals  has  already  been  worked  out. 


then  we  oan  find  that 


where  aAQ  Is  an  arbitrary  constant  of  Integration.  This 
Is  Just  relation  (3 .2) 

To  do  derivatives,  we  introduce  the  notation 

0 I J.  *«T«ry)J  - %o  a«  u ^ 


Then  by  relation  (3»l) 


By  equating  coefficients,  we  arrive  at 


where  a. 


Multiplication  can  be  done  using  the  relation 


aJ  fi 

If  <P  (y)  and  V'  (y)  » £0  Vn1^'  then 

2 Al 

- So  d*  Tnf/) 

_ £ £ o-i^;  />) 


„ 7 U^V  V L 

— i/=<> 

^ ^ Ft 

- •% 

*-o  r 0 


[71+/  (y)  + 71lyi  f/\} 


Truncating  this  sum,  it  can  be  written  as 

(86)  . _2  <khtt.sk  a±kc?+J-  jTfLktc 

(S.6)  dn  - ^ a.  + Cm  a + Cm  fa  2. 

Now  we  can  solve  equation  (7»52)  for  The 

left  side  is  Just  the  Orr-Sommerfeld  operator  with  °\ 
replaced  by  2°(.  The  equation  is  integrated  four  times, 
and  the  matrix  corresponding  to  the  left  side  can  be 
computed  using  the  same  subroutine  as  for  the  Orr-Sommerfeld 
equation. 


The  right  side  can  be  Integrated  explicitly  once  to 


obtain 


i.7j  „J>[  (^rr''t-  f'0  Dz<p'r'n] 


Using  our  subroutines,  the  corresponding  vector  can  be 
found  and  Integrated  three  times.  The  first  four  coeffi- 
cients are  set  to  zero,  since  these  equations  determine 
the  boundary  conditions.  The  result  is  a system  of 
inhomogeneous  linear  equations,  which  can  be  solved  to 
obtain 

Equation  (7»5l)  can  be  explicitly  integrated  twice. 


and  rewritten  as 


* 

(8.8,  p1  ^'"=  4 [->  ^ r* #n;n D *";"7 


The  vector  corresponding  to  the  right  side  can  be  determined 

jM  i) 

and  then  integrated  twice.  This  gives  us  F except 
for  the  first  four  coefficients,  which  are  determined 
from  the  boundary  conditions. 

Equation  (7.60)  can  now  be  solved  for  The 

definite  Integrals  can  be  evaluated  using  (3»1)  and  the 
formula 

(8.9)  T (+1)  * (tl)* 

I 

These  allow  us  to  compute  _ . . j -r  t1 

, \ ~r  / n i ^ r ^ J'H'tt'y ' Art)-}  //»*/  ( ft  / 

X,  J,  ^ \r«&>  - -^r, — J 

~ mzo  x L 

/>  T — ^2113  1 

— Cl/r\  (_  /*f^i  /H-t  J 

. ) n,  . _ 2^n_ 

(8.10)  J t JL^  QvT*  jlQo  + JfZx.  /XI  x-/ 

Combining  this  with  our  previous  subroutines,  and  ▲ 
can  be  explicitly  calculated. 


9,  NUMERICAL  RESULTS 

The  above  set  of  equations  was  solved  for  a variety 
of  values  of  and  R.  The  results  are  summarized  in 
Figure  2,  where  an  o represents  values  of  the  parameters 
where  a steady  state  solution  exists,  and  an  x values  for 
which  there  is  no  such  solution. 

Fig.  2. --Location  of  steady  state  solutions 


The  size  of  the  disturbances  associated  with  the 
steady  state  solutions  varies  widely.  Table  7 shows  the 
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tNtf1  'JJM*.  *1 


I i 


5*. 

amplitude  A of  the  steady  state  solutions  for  different 
values  of  and  R. 


TABLE  7 

VALUES  OF  R,°C  , AND  A 


R 

A 

R 

oC 

A 

: 100 

3 

1.57X10"1 

500 

13 

4.60xl0“6 

100 

4 

2 ,44xl0"2 

500 

14 

2 ,61x10-6 

100 

6' 

7.02xl0"3 

|t 

700 

9 

1. 76x10"* 

100 

7 

4.33x10  3 

✓ 

700 

12 

6.10x10“° 

100 

8 

3.20x10“-* 

/ 

700 

16 

1.18x10“° 

200 

5 

7.59x10-3 

800 

10 

4.83x10-3 

200 

6 

2,85x10"-* 

_4 

800 

12 

5.91x10“° 

200 

9 

3.51x10^ 

800 

15 

8.40x10”' 

200 

10 

2.19x10“* 

800 

17 

7.02x10*7 

200 

11 

4.12X10-4 

c 

1000 

15 

5.65X10”7 

o 

o 

9 

9.83x10  5 

, 

1000 

18 

2.45x10-7 

As  In  experiments,  the  smaller  the  Reynolds  number 
R,  the  larger  the  corresponding  amplitude.  The  smallest 
value  computed  was  R * 5 and  ^ « 2.  The  method  predicted 
a steady  state  solution,  but  with  an  amplitude  of  6.05* 
Since  our  formalism  assumed  small  disturbances,  this  does 
not  reflect  a realistic  solution. 


A 


It  was  not  feasible  to  determine  whether  there  exists 
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a largest  value  of  R for  which  there  is  a steady  state 
solution,  because  of  the  computational  expense.  To  obtain 
accurate  eigenfunctions  proved  more  difficult  than  to  obtain 
accurate  eigenvalues,  particularly  if  both  oi,  and  R were 
relatively  large. 

But  steady  state  solutions  are  possible  at  least  up 
to  R = 1000.  This  would  Imply  that  small  disturbances 
will  die  out.  However,  In  practice,  laminar  flow  can  be 
maintained  only  up  to  about  R = 750.  Since  in  an  actual 
apparatus  there  are  always  some  perturbations,  this  could 
be  due  to  the  fact  that  these  perturbations  can  no  longer 
be  reduced  below  the  transition  level.  Other  effects  not 
covered  by  our  formalism  could  also  be  at  work.  However, 
Just  considering  second  order  effects,  this  method  does 
correctly  predict  transition  to  turbulence  for  a large 
enough  perturbation,  which  the  linear  theory  does  not. 

Some  insicrhts  into  the  results  can  be  obtained  by 
considering  the  shape  of  the  disturbances.  Using  a second 
order  approximation 


(9.1) 


f (A.y,  ») 


7 £ +r  $ -f  f j. 


(9.2)  f (A ,y.&) 

-I (fn  t i»+]  + 3.  A 1 i <p’ fc>- * ■+■ 

Re«l  (<P(Vin)  L.oil-0-  -lr*y(4r*'x)  sinl*] 
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( (p  1 is  a real  valued  function*)  Then  since  u = ^ 2^^ 

and  ▼ * * 

(9.3)  »(A.y.*>  - y * ZaZ**'1  *"*>  C”  + 

. <Pr,:v  ’>«+: +^it>* //;j’ 

+ t>(ltc«l  0fX'13)c.»  i'?  - l>  (l«‘)0n'V)  >;n 

(9.4)  Y(A.y.«->  - -XA[(fteJ  (S f>ia)*r*+ -/$*»?  (t"'V <■>**] 

(Heal  0 >V,U)  S/o  1-6  -tJJnHf  (pn'^)  CO*  x 

If  we  ignore  the  A2  terms,  we  have  the  Orr-Somraerfeld 
linear  solution* 

Choosing  » 0,  this  simplifies  to 

2 A I)  J aO,))  ) 

(9.5)  u(A.y,0  ) * y + A OiKeair  y 

+ £ A*  Z t>0";"  + t>  (A**  / <P '*  *V ] 

(9.6)  y(A,y,  0 ) - -2/9 Imy  0"'0-  0 11 

Using  the  definition  Tn(cos  T ) * cos  , these 
vectors  can  be  computed.  This  Is  done  most  easily  by 
choosing  values  of  % between  0 and  T[  * Then  y ■ cos  If 
and  Tn(y)  « cos  n f , 

The  resulting  solutions  are  graphed  for  the  case 
R * 200  and  6<  ■ 5*  Figure  3 gives  u and  v for  the  linear 
case,  and  Figure  4 gives  the  non-linear  solutions* 

For  the  linear  case,  the  disturbance  is  concentrated 


near  the  boundary  y ■ 1.  There  is  of  course  a symmetric 
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solution  concentrated  near  the  boundary  y = -1,  corresponding 
to  an  eigenvalue  with  negative  real  part.  For  larger  R, 
the  disturbance  is  even  closer  to  the  boundary,  but  retains 
the  same  shape. 

In  our  case,  we  are  multiplying  by  an  amplitude 
factor  A.  As  we  have  seen,  this  is  related  to  the  amount 
of  change  in  the  fluid  velocities  with  resnect  to  x and  t. 

The  size  of  A is  immaterial  if  the  disturbance  is  small 
enough  for  the  linear  theory  to  hold. 

The  actual  magnitude  of  the  fluid  velocities  depends 
on  pn)0 

as  well  as  on  A.  We  have  used  the  normalization 
(ft  'l>/0)  = 1,  but  the  size  of  this  eigenfunction  still  varies 
greatly.  If  ck  is  fixed  and  R Increases,  the  size  of  this 
function  Increases  slowly.  But  if  R is  kept  fixed  and  c< 
Increases,  rapidly  becomes  larger.  For  Instance, 

if  R = 200  and  A is  kept  fixed,  the  maximum  value  of  u 
for  of  = 12  is  almost  a thousand  times  as  large  as  the 
maximum  value  of  u.  for  = 5.  So  for  large  °<  , A must 
be  chosen  quite  small  if  the  linear  theory  is  to  apply. 

The  shape  of  the  non-linear  solution  is  completely 
different  from  the  linear  solution.  However,  this  same 
basic  shape  is  ag^ln  retained  for  different  values  of 
and  fl.  In  fact,  in  most  cases,  the  magnitudes  are  also 
quite  similar.  As  we  have  seen,  for  a given  fixed  R, 

A decreases  rapidly  as  increases.  However,  this  does 
not  decrease  the  difference  between  the  non-linear  solution 
and  the  laminar  solution,  since  the  eigenfunctions  are 
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getting  large  rapidly. 


The  value  of  A la  now  important,  since  the  given 


second  order  solution  only  exists  for  the  appropriate 


amplitude.  When  we  multiply  by  A,  as  we  have  been  doing. 


the  linear  solution  becomes  closer  and  closer  to  the 


laminar  solution  as  R increases.  But  the  second  order 


solutions  remain  distinct  from  the  laminar  case. 


The  transition  from  the  linear  to  the  second  order 


solution  cannot  be  obtained  by  our  method.  However,  we 


know  that  for  small  disturbances  the  linear  theory  holds. 


Somewhere  well  before  the  transition  amplitude  non-linear 


effects  would  begin  to  occur.  For  these  intermediate 


disturbances  tho  shape  would  begin  to  change  from  the  linear 


to  the  non-linear  form.  By  the  time  the  transition  ampli- 


tude is  reached,  the  non-linear  effects  have  drastically 


changed  the  shape  and  the  size  of  the  solution.  So  our 


theory  shows  a large  difference  between  the  linear  and 


non-linear  solution,  which  explains  how  the  transition 


to  turbulence  can  occur. 


LIST  OF  REFERENCES 


BETCHOV,  R.  and  CRIMINALE,  W.  196?.  “Stability  of  P&rallel 
Flows.1*  Academic  Press  Inc. 

DAVEY  A.  1973.  On  the  stability  of  plane  Couette  flow  to 
infinitesimal  disturbances.  J.  Fluid  Mech.  vol.  57 » 
pp.  369-80. 

FORSYTHE,  G.  and  MOLER,  C.  1967.  "Computer  Solution  of 
Linear  Algebraic  Systems."  Prentice-Hall,  Inc. 

FOX,  L.  and  PARKER,  I.  1972.  "Chebyshev  Polynomials  in 
Numerical  Analysis."  Oxford  University  Press 

GALLAGHER,  A.  and  MERCER,  A.  1962.  On  the  behaviour  of 
small  disturbances  in  plane  Couette  flow.  J.  Fluid  Mech. 
vol.  13,  pp.  91-99. 

GARY,  J.  and  HELGASON , R.  1970.  A matrix  method  for 
ordinary  differential  eigenvalue  problems.  J.  Comp.  Phys. 
vol.  5*  PP»  169-87. 


ORSZAG,  S.  1971.  Accurate  solution  of  the  Orr-Sommerfeld 
stability  equation.  J.  Fluid  Mech.  vol.  50,  pp.  689-703* 

OSTROW3KI,  A.  1958-59.  On  the  convergence  of  the  Rayleigh 
quotient  iteration  fot  the  computation  of  characteristic 
roots  and  vectors.  Arch.  Rat.  Me;h.  Anal.  vol.  1,  pp. 
233-^1.  vol.  2,  pp.  423-28,  vol.  3,  pp.  325-40,  vol.  3. 
pp.  341-47,  vol.  3,  pp.  472-81,  vol.  4,  pp.  153-65. 

REICHARDT,  H.  1956.  Uber  die  Geschwindlgkeitsverteilung 
in  elner  geradlinigen  turbulenten  Couettestromung.  Z. 
Angew.  Math.  Mech.  vol.  36,  supplement,  pp.  26-29. 


REYNOLDS,  W.  and  POTTER,  M.  1967.  Finite-amplitude  insta- 
bility of  parallel  shear  flows.  J.  Fluid  Mech.  vol.  27, 
pp.  465-92. 

WILKINSON,  J.  1965.  "The  Algebraic  Eigenvalue  Problem." 
Clarendon  Press -Oxford 


WILKINSON,  J.  and  REINSCH,  C.  1971.  "Linear  Algebra." 
Sprlnger-Verlag 


61. 


